Tracing the sources of cosmic rays with molecular ions 
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ABSTRACT 

The rate of ionization by cosmic rays in interstellar gas directly associated with 7- 
ray emitting supernova remnants is for the first time calculated to be several orders of 
magnitude larger than the Galactic average. Analysis of ionization-induced chemistry 
yields the first quantitative prediction of the astrophysical emission line spectrum, 
which should be detectable together with lines. The predicted coincident observation 
of those emission lines and 7-rays will help prove that supernova remnants are sources 
of cosmic rays. 

Subject headings: supernovae: general — cosmic rays — acceleration of particles — 
stars: winds, outflows — shock waves — radio continuum: general 



1. Introduction 

One of the central questions in modern physics concerns the origin of cosmic rays (CRs). The 
CRs themselves do not directly reveal their sources as they interact with cosmic magnetic fields 
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and particles on their way through the Universe (jStrong et alj 120071 ). One method of identifying 
CR sources indirectly is the search for high-energy photon or neutrino signals. When protons of 
energy E > 200 MeV interact with photon or matter fields in the vicinity of the CR production 
site, charged and neutral pions are produced, which then emit high-energy neutrinos and photons. 
Those neutral particles can be used to trace CRs, as they point back to the sources. This approach 
is liin ited by (a) the low interaction cross sections that make neutrinos difficult to detect (iBecker 
20081'). and (h) the competition with processes like inverse Compton scattering and bremsstrahlung 
(|Schlickeiseiil2002l ) . It is therefore crucial to use multiwavelength approaches to identify CR sources 
through correlation studies, and not to rely solely on the observation of high-energy photons. 

The i llumination of rnolecula r clouds by CRs fr om supernova remnants (SNRs) was first sug- 
gested by iBlack and Fazid (jl973l ): iMontmerld ()1979l ). Most recently, observations performed with 
the Fermi Telescope have revealed high-energy emission from several regions associated with a SNR 
located in close proximity to molecular clouds (MCs). In parti cular, en iission at 0.2 to 20 GeV 



has been detected towards t he SNRs W51C f Abdo et al.ll2009lV W28 f Abdo et al.l l2010al'l. W44 



(|Abdo et al.ll2010bl ). W49B (|Abdo et al.ll2010cl l and IC443 (|Abdo et al.ll2010dl). SNRs are one o f 
the primary candidates for CR acceleration, as discussed for example bv iBiermann et al.l (j2009l ). 
The observed gamma-ray signatures could therefore be interpreted as vr''— decay products from 
high-energy CR interactions with the MC. The data f or these particular objects actually seem to 
favor a hadronic scenario over a leptonic one, see e.g. (jAbdo et al.ll2010bl ). Due to the ambiguity 
of the signal, however, other tracers of CRs are needed to have a final proof. 

In this paper, we quantify a new method to identify SNRs as CR sources, based on CR-induced 
ionization. When CRs are accelerated in a SNR and then injected into a nearby MC, the ionization 
rate in the ambient gas is expected to be enhanced compared with the ionization caused by the 
general background of Galactic CRs. The immediate products are hydrogen ions, H"*" and Hg", 
which react rapidly with II2 and oxygen to form molecular ions. Ions like Hj^, OH"'", 1120"*", and 
HsO"*" are spectroscopically observable at infrared and s ubmm wavel e ngths and their a bundances 
are predictable as functions of the ionization rate (e.g. (iBlackl Il998l : IMontmerld l2009l )). Among 
the 18 molecular cations now identified in the interstellar medium, it is the reactive ions Hj", 



OH"*" , and H^O"*" that are expected to be the most specific tracers of the CR ionization rate (jBlack 



2007l\. and it was only during the last year that the latter two spec i es have been detected 
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of 011"*" and 1120"*" absorption in diffuse interstellar gas of low molecular fraction (H2/II < 10%) 
serves to calibrate the backgro und ionizing frequency of interstellar h ydrogen at a level Ccai = 
(0.6 - 2.4) • lO^iSg-i pgr atom (jCerin et ahlboiol : iNeufeld et ahlboiOal M. 
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2. Ionization in the direct vicinity of supernova remnants 

Here, we estimate the enhanced interstellar ionization rates expected from the observed fluxes 
of high-energy gamma rays toward SNRs and predict the abundances of reactive molecular ions 
to be expected in associated MCs. We show how correlations between GeV-TeV emission and 
vibrational and rotational spectra of molecular ions could be used to identify the sites of CR 
acceleration. Observational tests will be possible through astronomical spectroscopy at infrared 
and mm/submm wavelengths. CRs of kinetic energies between E'j, ~ 1 MeV and ~ 1 GeV ionize 
the interstellar medium, while those at lower energy hardly escape their sources. We investigate 
five SNRs located in or close to MCs, which are detected in gamma-rays. 



2.1. Primary particle spectra 

The CR spectra, further used to calculate the ionization rate, can be deduced from the mea- 
sured gamma-ray flux at Earth: 



dE^ dt dAEarth 

Assuming that the photons arise from hadronic interactions between the SNR and the MC, a 
primary CR spectrum can be derived. We assume a smoothed broken power-law for the primaries: 

dNp _ / p y f f p \ Y^' 

dEp dt dAsource " 1 GcV C"! / V ^ V W J ' 

where Op is the proton spectrum normalization and pbr the location of the spectral break. For all 
sources, this shape of the primary spectrum reproduces the sha pe of the observed photon spectrum. 
Our results generally follow the original calculations by Fermi, ( Abdo et al. 2009 . 2010al lbllaldl) with 



slight modifications of the smoothing of the power-law, only giving minor changes to the result. 



We use the parametrization of the CR interaction model as presented in lKamae et al.l ()2006l ) to 
derive the primary spectrum from the gamma-ray observation. We use a fixed value of the hydrogen 
density, nn = 100 cm~^ and determine the CR energy Wp needed to provide the observed gamma- 
ray flux. The volume of the emitting region V enters the normalization of the proton spectrum. It 
can be derived from gamma-observations and is taken from the Fermi publications. The spectral 
parameters s, As and phr, the total energy budget of interacting protons, Wp and the volume V 
are listed in Table [TJ 



2.2. Ionization 



To calculate the ionization rate, we use the parametrization from iPadovani et al.l (|2009l ). in- 
cluding ionization by CR electrons, CR protons and electron capture by CR protons. Only the 
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W51C 


W44 


W28 


IC443 


W49B 


s 


1.5 


1.74 


1.7 


2.09 


2.0 


As 


1.4 


1.96 


1.0 


0.78 


1.4 


Pbr [GeV c"i] 


15 


9 


2 


69 


4 


Wp [10^9 erg] 


9.0 


12.1 


3.3 


28.0 


35.0 


V [cm'^] 


3.3 • 10^° 


4.2 • 10^9 


3.2 • 10^^ 


4.2 • 10^9 


6.3 • W^^ 



Table 1: Parameters used to derive the primary cosmic ray spectrum from the observed photon 
spectrum . The total energy budget of the sources and the spectral shape are determined from the 
observed gamma-ray spectrum, the interacting volume is deduced from the emitting region (see the 
original Fermi papers). 



direct ionization by primary protons turns out to be significant, the other contributions being at 
least one order of magnitude smaller. The ionization rate of H2 by primary protons reduces to 

C^^=/ -^<{Ep)dEp, (3) 

with dNp/dEp as the CR spectrum at a given kinetic energy Ep. The upper integration limit is set 
to i?max = 1 GeV, as the rapidly decreasing cross section makes contributions at higher energies 
negligible. 

There are two main factors affecting the ionization rate: 



(1) The CR spectrum at Earth as measured below ~ 1 GeV is modified from the spectrum at 
the source by transport effects of the ISM and the solar environment. The particle spectrum 
at the source can be derived from gamma-ray data up to approximately 1 GeV, but it must 
be estimated theoretically at lower energies. The simplest approach is to extrapolate the 
spectrum towards lower energies, assuming that the acceleration mechanism works the same 
way at MeV energies. However, low-energy partic les might follow a different acceleration 
scheme and produce a flatter or curved spectrum, see lBlasi et al.l (|2005l ) and references therein. 
This might lead to an overestimation of the spectrum at the lowest energies. However, it is 
not clear yet at what energy this will start exactly. 

(2) The choice of lower energy threshold -Emin is related to propagation effects. The range of a 1 
MeV proton in hydrogen is 8.5x 10~^ g cm~^, which corresponds to a projected column density 
of 5.1 X 10^*^ H cm~^. This is easily able to sup port detectable abundances of 0H+ and H2O"'" 
in the diffuse molecular gas of the Galaxy. In llndriolo et al.l ( 20091 ). a threshold of -Bmin = 2 
MeV is recommended for diffuse clouds with hydrogen number densities uh < 10'^ cm"'^, and 
Emin = 10 MeV for dense clouds, where nn > 10'^ cm^^. Close to a CR accelerator, at an 
average density of 100 cm~'^, values above 2 MeV are reasonable choices to account for the 
stopping range. 
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We chose the fohowing approach to account for these uncertainties: we first assume that the 
spectrum continues as a power-law down to MeV energies. We then introduce an artificial and 
sharp cutoff in the spectrum at a minimum kinetic energy E'mim below which the CR flux is set 
to zero. We then vary the minimum energy between Smin = 1 MeV and 100 MeV. As the CR 
ionization cross section decreases with energy, increasing the threshold will remove much of the 
potential signal, having an even more extreme effect than flattening the CR spectrum. Increasing 
the threshold at the same time accounts for considering different stopping ranges. At this stage, 
due to the lack of knowledge of the lower part of the spectrum, we will calculate the possible range 
of ionization intensity rather than one specific value. 



— W49B — 


W28 


— -[C443 — 


W51C 


— «44 


stopping range 




E„,„ [eV] 

Fig. 1. — Dependence of the ionization rate, normalized to the average Galactic value of Ccai '■- 
2 ■ 10-^^ s-^ on £;min. 



Figured] shows the ionization rate, normalized to Galactic average, versus the minimum energy. 
For comparison, the stopping range is shown as a function of -Emin- For the most optimistic scenario, 
still compatible with current knowledge, E'min = 2 MeV, ionization rates range from 110 to 4 • 10^ 
times the Galactic average. An energy cutoff at 100 MeV leaves W49B as an outstanding source 
with a factor o f 4 • 10^ above Gal actic average. The central reason is the large CR energy density 



in this source (lAbdo et al.ll2010d ). At i^min = 100 MeV, only W51C drops to values below Ccai- 
All other sources remain a factor of a few to ten above the background. 

Note that we do not include the additional ionization by energetic secondary electrons being 
produced in the primary ionizing event. These secondary ionizations increase the total rate by 
factors of a few (cf. Pa dovani et al. 2009), bu t the details are sensitive to the relative abundances 
of H, H2, e~, and He (jPalgarno et al.Ml999l ). The escape of CRs from their sources and their 
propagation into su rrounding material a r e also affected in complicated ways by interaction with 
magnetic fields (see iPadovani and Gallil (|201lh for a recent discussion for dense molecular gas). 
These issues will be addressed in future work. 
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3. Cosmic ray induced molecular ions 

As shown below, the observable line intensities of molecular ions spawned by CR ionization 
depend upon the ionization rate integrated through a column of molecular gas 



[ C^^n(R2)di ^ C"2^(H2)L s~^ cm'^ 
Jo 



which we call the cosmic ray exposure. The number density of interstellar H2 is n(H2) and the 
integration is over path length ^ = to L. Because the CRs of shortest range die off at small i, 
protons at energies significantly greater than E'mm contribute to the exposure. In principle, the 
exposure should be computed through an explicit depth- integration of the discrete energy loss. The 
SNRs discussed in this context have a CR exposure of the order of ?] ~ 109 S-1 cm ^ or higher 
(see Fig. [1]), which is > 10^ times larger than the exposures that characterize diffuse clouds in the 
Galactic disk. In general, ionization by X-rays can produce similar chemical signatures. We expect 
the X-rays from SNRs to be much less important than CRs because the X-ray luminosities are 
relatively small and their depths of penetration are too short to yield high values of the exposure. 



The outlines of interstellar ion chemistry have been understood for almost 40 years (jHerbst and Klemperer 



I973I ) and the rate coefficients of all the crucial reactions are well determined. In steady state be- 
tween the rates of source and sink reactions, the number densities of the transient ions O'^, OH"*", 
II2O"'", and HeH"'" are expected to be proportional to (^h- In fully molecular regions, n(H) << n(H2), 
of low ionization x{e~) < 10^^, competing processes rarely interrupt the sequence of H-atom ab- 
straction reactions so that almost every ionization of hydrogen leads to formation of , Hj^ , OH"'" , 
112©"*", and HaO"*". The time-scales to achieve steady state are very short, ~ l(lOOO cm~'^/n(H2)) 
year. 

The abundances of interstellar molecules can be computed from a system of rate equations 
for various values of the physical conditions (density, temperature, ionization rate) and gas-phase 
elemental abundances, either in the steady-state limit or as a time-dependent solution. However, 
the abundances alone do not characterize observable signatures, because the distributions of the 
molecules among their quantum states tend to be far out of equilibrium. The departures from 
equilibrium will be most extreme for the most reactive molecular ions, because elastic and inelastic 
collisions that would otherwise thermalize them are tvpic a.lly no t faster than the reactive collisions 



with the main collision partner, H2, that destroy them (iBlackl llQOS ). We illustrate the coupled 
effects of ion chemistry and non-equilibrium excitation with results from a reference model of 
molecular gas at density n(H2) = 100 cm~^ exposed to an ionization rate = 10^^^ s~^ over a 
path length L = 10^^ cm = 3.24 pc. The excitation and transport of spectrum lines are computed 
through use of the Radex program, which can incorporate formation and destruction rates state- 
by-state for reactive mo lecules along with all relevant radiative and inelastic collision processes 



(|van der Tak et al.l 120071 ). The kinetic temperature of the molecular gas is assumed to be 150 
K, somewhat elevated owing to enhanced heating by cosmic rays. All line widths are taken to 
be 9.394 km s~^. Full details of these calculations will be presented elsewhere. Because H^ has 
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no permanent dipole moment, its vibration-rotation and pure rotational spectrum occurs only as 
weak electric quadrupole transitions. Upon birth by ionization of H2, is endowed with high 
vibrational excitation but it retains the relatively cold rotational distribution of its parent. At 
gas densities below 10^ cm~'^, excited decays to the ground state by radiation more rapidly 
than it is destroyed or excited by collisions. As a result, the emitted spectrum in our reference 
model is directly predictable from the val ue of the CR exp osure, rj = 10^ s~^ cm~^, and the 
known quadrupole transition probabilities (jPosen et al.Ml983l ). Owing to symmetry-breaking in 
the last few bound vibration-rotati on levels near the dissociation limit, Hj" exhibits a few dipole- 
allowed radio- frequency trans i tions (Brown and Carrington 2003 : Howells and Kennedvl[l991 : Moss 



1993 



Bunker and Moss 



200c 



Critchlev et al 



2001 



These transitions are excited in our model by 
the small fraction of ionizations that put 112" into levels {v = 19, J = 0, 1) of the Iscjg ground 
electronic state. Transitions t o three bar ely bound levels {v = 0, J = 0,1,2) of the 2pau state 
occur with high probabilities (jMosslll993l ) that partly compensate for the low populations of the 
initial states: we predict that these transitions at 17.6, 52.9, and 96.4 GHz will appear as weak 
masers in regions of high CR exposure. The electronically excited levels also radiate rapidly to 
lsag{v = 17 — 18, J = — 3) in transitions at mm/submm wavelengths. Although these symmetry- 
breaking transitions will probably be undetectably weak in sources like those considered here, they 
may become observable in the more extreme environments in active galactic nuclei. 

The high excitation of is partly transferred to Hj^ in the next chemical step and Hj]" can also 
be excited by collisio ns several times during its lifetime. The distortion-induced rotatio nal transitio n 
(jPan and OkJl98d ) of H3 (J, K) = (4, 4) (3, 1) easily suffers population inversion (lBlacklll998l ). 
a necessary condition for a maser. The transition frequency is poorly known (approximately 217.7 
GHz), but the line is predicted in our reference model to appear as a weak maser with a peak 
Rayleigh-Jeans brightness temperature of 10 mK, an intensity that is detectable with existing radio 
telescopes. The intrinsically weak, far-infrared rotational transitions of H^^ are excited mainly by 
collisions at the kinetic temperature of the gas, while the infrared vibration-rotation transitions of 
both H^ and are excited in the formation process with superthermal intensity distributions. In 
general, the rotational line intensities of Hj^ will depend upon the local density and temperature; 
however, the entire spectrum of H^ is rather insensitive to these physical conditions. We calculate 
explicitly the formation rate of HeH+ by vibrationally excited H^ and predict the spectrum that 
results from its initial, superthermal excitation. The predicted spectra are illustrated in Figure [2j 
At infrared wavelengths the strongest lines of Hj^ appear between 3 and 4.6 ^m while Hg" lines 
of comparable intensity appear between 4 and 6 ^m. The strongest vibration-rotation lines of 
HeH+ in this region are several hundred times weaker. In the far-infrared region, the thermally 
excited lines of Hj^ are quite strong while the lines of HeH+ are again considerably weaker. In 
the far-infrared, thermal emission by interstellar dust particles will produce continuous radiation 
with typical surface brightness > 0.1 Jy nsr~^ at a wavelength of 100 ^m; even so, it appears 
that the strongest rotational lines of might be observable if the spectroscopic resolving power 
is at least 100 (resolution of 1 /.tm at 100 /.tm wavelength). For comparison, the weak Hg" feature 
at 30.725 /.im has an intensity approximately 1/2 as large as that of the H2 J = 2 — )■ line at 



0.02 



0.015 



0) 0.01 



5x10 



CN(H2) = 10^ at n(H2) = 10^ 



jullJUjiiLil 
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200 



Fig. 2. — Predicted emission spectra of (red in online version) and Hj]' (black online) for the 
reference model. The surface brightness is displayed in Jy nsr~^ where 1 Jy = 10~^^ W m~^ Hz"^ 
and 1 nsr = 10^^ steradian of solid angle. In the figure, rj = 10^ s"^ cm"^ is used as a typical value 
in this context (see text). 



28.226 fxm. We have also calculated the spectra of the OH^ ions, n = 1 — 3. These ions form 
mainly as a result of reactions of thermalized Hjj" or H+ with oxygen and subsequent reactions 
with H2. The excess energy (enthalpies of reactants minus products) of each reaction will go partly 
into vibrational and rotational excitation of the product, so that the OH^ ions will also exhibit 
superthermal excitation. A predicted signal of enhanced production of oxygen-bearing ions will 
be strong absorption in rotation-inversion transitions arising in highly metastable states of HaO"*" 
at frequencies 1.6 to 1.9 THz. Absorption spectroscopy of interstellar Hj^ has alread y suggested 



enhanced rates of ionization throughout the central molecular zone o f the Galaxv (jGoto et al 



2008; 



Geballe and Qkall20ld ) and in two locations near the SNR IC 443 (jindriolo et al.ll2010al ). We 



predict that emission lines of both and Hg" will be detectable and that their intensities can be 
used to measure the CR exposure. 
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4. Summary 

We have calculated the associated cosmic ray ionization rate in interstellar matter near SNRs. 
While the exact value will depend on the exact density and spectral shape, we expect ionization 
to be significantly above Galactic average. The cosmic ray exposure of the interstellar hydrogen 
(ionization rate times gas density integrated over path length) for the SNRs W51C, W28, W44, 
IC443 and W49B is expected to be of the order ofrj = 10^ cm~^. This permits the construction 
of a coupled model of the ion chemistry and the excitation of reactive molecular ions that might 
serve as specific tracers of cosmic rays near their sources. Predictions are made of the intensities 
of rotational and vibration-rotation transitions of the reactive ions including the first quantitative 
prediction of an astronomical spectrum of H^. It appears that superthermally excited emission lines 
of and Hj^ may be detectable at wavelengths between 3 and 6 /im. Strong thermal emission 
of Hg" may be measurable in the far infrared at 50 to 150 fim wavelength. W49B seems to be 
an outstanding candidate for the detection, with predicted ionization rate more than 4 orders of 
magnitude above Galactic average. Future measurements by e.g. Herschel and ALMA provide the 
opportunity to make observations of ionization signatures in SNR/MC systems and analyze their 
correlation to the gamma-ray emissior i. A recent observati on suggests an enhanced excitation and 



abundance of ammonia toward W28 (iNicholas et al.ll201ll ). Using Radex, it is possible to explain 
the observed spectrum as CR induced with an ionization rate of 4.5 ■ 10~^^ s"^, so a factor of > 10^ 
above Galactic average. Further, Hg" emission observed from the direction of IC44 3 indicates an 
enhanced ionization level of 2 ■ 10~^^ s~^, so a factor of 10 above Galactic average ( Indriolo et al 



201 Obi ). In the future, more observations need to be performed in order to explicitly search for the 
expected signatures of and Hj^. Treating the propagation of the CRs in detail will allow us to 
give even more precise predictions of the expected emission regions. 
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